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ABSTRACT 

We have succeeded in obtaining exact configurations of irrotational binary systems 
for compressible (polytropic) equations of state. Our models correspond to binaries 
of equal mass neutron star systems in the viscosity free limit. By using the obtained 
sequences of stationary states, the evolution of binary systems of irrotational neutron 
stars due to gravitational wave emission has been examined. For inviscid binary sys- 
tems the spin angular velocity of each component in a detached phase is smaller than 
the orbital angular velocity at a contact phase. Thus the irrotational approximation 
during evolution of binary neutron stars due to gravitational wave emission can be 
justified. Our computational results show that the binary will never reach a dynam- 
ically unstable state before a contact phase even for rather stiff polytropes with the 
index N <; 0.7, as the separation of two components decreases due to gravitational 
wave emission. This conclusion is quantitatively different from that of Lai, Rasio & 
Shapiro who employed approximate solutions for polytropic binary systems. 
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1 INTRODUCTION AND SUMMARY 

Coalescing stages of binary neutron star systems due to 
gravitational wave emission have been considered to be 
one of the most important targets of the advanced gravi- 
tational wave detectors (LIGO /VIRGO /GEO /TAMA, see 
e.g. Abramovici et al. 1992 and Thorne 1994). Observations 
of final stages of binary neutron star systems will provide 
us a large amount of new information about macroscopic 
quantities such as the mass and the spin of neutron stars as 
well as about microscopic characters such as the viscosity 
and the equation of state of neutron star matter ( Cutler et 



Co(x) = V x v(x) 



(1) 



al. 1995) 



Kochanek ( |l992j ) and Bildsten & Cutler fll992j ) pointed 
out that viscosity in neutron stars may not be so effective as 
to synchronize the spin and the orbital angular velocity on 
a time scale of evolution due to gravitational wave emission. 
For such inviscid fluids, Ertel's theorem holds: the ratio of 
the vorticity vector Co in the inertial frame to the density 
p of a fluid element, Co/ Pi IS conserved even under the exis- 
tence of a potential fo rce such as the gravitational radiation 
reaction (Miller 1974). The vorticity vector in the inertial 
frame is defined as 



where x and v(x) are the position vector of the fluid element 
and the velocity field seen from the inertial frame, respec- 
tively. Since the spin angular velocity of each component of 
a detached binary system is expected to be at least 50 times 
smaller than the orbital angular velocity at an almost con- 
tact state of the binary system (see the references above), 
we may consider that real close binaries have irrotational 
configurations. In other words, the vorticity of each compo- 
nent star seen from the inertial frame of reference can be 
neglected, i.e. Co = 0. 

It is not, however, an easy task to construct consistent 
models of stationary configurations of compressible stars 
such as binary systems with arbitrary spins (see e.g. Uryu & 
Eriguchi (1996) and references therein). Therefore compress- 
ible binary configurations in equilibrium states have been 
investigated only for synchronized binary systems in Newto- 
nian gravity (e.g. Hachisu & Eriguchi 1984; Hachisu 1986) or 
in post-Newtonian gravity (Shibata 1994, 1997) up to highly 
deformed configurations. As for configurations whose spins 
are different from the orbital angular velocities, there are 
only approximate solutions by Lai, Rasio & Shapiro (1993a, 
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1994a (hereafter LRS1), 1994b (LRS2)). They employed tri- 
axial ellipsoidal polytropes for deformed binary states in 
Newtonian gravity and discussed evolutions of binary neu- 
tron stars. Since real configurations are no more ellipsoidal 
when the binary comes near to a contact phase, their results 
cannot give quantitatively correct values even in Newtonian 
gravity. 

Recently we have formulated a scheme to compute irro- 
tational binary configurations composed of two compressible 
stars with equal mass and developed a numerical code to 
solve exact conf igurations for such binar y systems in New- 
tonian gravity (Uryu & Eriguchi 1997b). Using this newly 



developed code, we have constructed stationary sequences 
of irrotational binary systems. Obtained models have been 
used to investigate evolutions of binary systems due to gravi- 
tational wave emission. Our computational results show that 
dynamical instability will not set in before a contact phase 
for polytropes with the polytropic index TV <; 0.7 as the sep- 
aration of two components decreases. For polytropes with 
the index in this range, dynamic al instability will occur on 
a so-called 'dumbbell' sequence (Eriguchi & Hachisu 1981: ) 
for irrotational self-gravitating fluids. This result is different 
from that of Lai, Rasio & Shapiro (1993a, 1994a, 1994b) who 
concluded that dynamical instability sets in on polytropic 
binary sequences with N ^ 1.2 by using ellipsoidal approx- 
imations of deformed configurations. In this Letter we will 
briefly summarize our new result and discuss its physical 
relevance in evolution of binary neutron star systems. 



2 FORMULATION OF THE PROBLEM 

We consider stationary structures of polytropic binary stars 
without viscosity in the rotating frame. We assume that the 
binary is composed of equal mass components whose vortici- 
ties equal to zero in the inertial frame of reference. For such 
irrotational cases, we can introduce the velocity potential 
$(x) as follows: 



(2) 



where v is the velocity vector in the inertial frame. The Eu- 
ler equation of fluids can be integrated to the generalized 
Bernoulli's equation as follows in the rotating frame of ref- 
erence: 



<9$ 1 
_- (n xr).V* + - 



v$ r + 



dp 
P 



+ $ = /(t).(3) 



where p, cj>, fl and r are the pressure, the gravitational poten- 
tial, the orbital angular velocity vector of the binary which is 
identical to the angular velocity vector of the rotating frame 
relative to the inertial frame and the position vector of the 
fluid element of the star from the rotational center, respec- 
tively, and f(t) is an arbitrary function of time. We note 
that the steady velocity u in the rotating frame of reference 
is related to the velocity v in the inertial frame as follows: 



CI x r . 



(4) 



In actual computations we use the polytropic relation 
and the Emden function defined as: 



Kp 1+1/N 



KQ 



N+l 



(5) 



where O is the Emden function and A' is a constant. Since 
the configuration of the binary is assumed to be stationary 
in the rotating frame, we can set 



<9$ 

— — = and f(t) = C — constant, 
dt J w 

Then equation ^ becomes as 

1 



e 



K(N + 1) 



(6) 



■<f> + C].(7) 



The equation of continuity is expressed by using the 
velocity potential as follows also in the rotating frame: 

V 2 $ = N(fl x r - V$) ■ (8) 

where the assumption of a stationary state in this frame has 
been taken into account as follows: 



dp 
dt 



o. 



(9) 



After substituting the following integral expression for 
the Newtonian gravitational potential 



= -G 



-d r 



(10) 



into equation (m), where integration is performed over the 
whole stellar interior V , there remain two unknown physical 
variables in the basic equations, i.e. the Emden function O 
and the velocity potential $. We can regard equation (m) 
as the equation for the variable O and equation (^) as that 
for the variable $. The boundary conditions for these two 
variables are as follows: 



(V$ - f2 x r) ■ n = 0, 

6 = 0, 



on the stellar surface, (11) 
on the stellar surface. (12) 



These basic equations are transformed into the surface fit- 
ted coordinate system (Uryu & Eriguchi 1996) and solved 
iteratively by using the self-consistent field method (Hachisu 
1986). The detailed nu merical method will be e xplained in 
the forthcoming paper ( Uryu fc Eriguchi 1997b ). 



3 EVOLUTIONARY SEQUENCES 

We have computed several sequences along each of which 
the mass, the entropy (K) and the circulation (zero) are 
kept constant. Thus these sequences can be regarded as the 
quasi-statically evolving sequences of inspiraling binary neu- 
tron stars due to gravitational wave emission. In figures 1, 2 
and 3 we compare our results with those of the ellipsoidal ap- 
proximation of the irrotational compressible binary systems 
of LRS1. The physical quantities of stationary sequences of 
irrotational equal mass binary systems with TV = 0, 0.5 and 
1 are shown in figures 1, 2 and 3, respectively. In these fig- 
ures, terminal points at smaller separations or smaller an- 
gular momenta correspond to models at contact stages of 
binary stars. The normalized values for the total angular 
momentum J, the total energy E, the angular velocity Q 
and a separation between the mass centers of two compo- 
nent stars da are defined as follows, 



J = 



J 



(GMZRo) 1 / 2 



E 



E 



GM 2 /R 
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Figure 1. Physical quantities of stationary sequences of irrota- 
tional binary states for N = incompressible fluid stars, (a) Total 
angular momentum as a function of a binary separation, (b) Total 
energy as a function of a binary separation, (c) Orbital angular 
velocity as a function of the total angular momentum. Dashed 
and solid curves show the results of LRS1 and our present results 
of irrotational binary stars, respectively. Dots show the results 
computed by using different numerical method (Uryfl & Eriguchi 



1997a). pee text about the normalization factors tor each quan- 
tity. 



and cLg 



Ro 



(13) 



where G is the gravitational constant, M is the mass of 
one component star, Ro is the radius of the spherical star 
with the same mass M and the same polytropic index N. 
The quantity po is defined as po = M/(AtvRq/3). This 
normalization is the same as that of TRST Results are in 
good agreement with each other for models with larger sep- 



Figure 2. Same as figure 1 but for the results of LRS1 and our 
present results with the N = 0.5 equation of state. 



arations. Differences at the separation ^ 3.5 are at worst 
& 0.5% for any quantity. In figure 1 we also show the re- 
sults computed by our new computation al method for in- 



compr essible binary with internal flow (Uryu & Eriguchi 



1997a). Relative error between them are less than ~ 0.5% 



everywhere. These figures ensure the accuracy of our present 
results since the two independent computational methods 
give the same results for N = case even near the contact 
phase where the stars are significantly deformed by tidal 
field. 

On the other hand, for smaller separations, differences 
between LRS and our results become evident. For ellip- 
soidal approximations, curves with N & 1.2 have turning 
points where the total angular momentum and the total en- 
ergy attain their minimum values. In our irrotational binary 
models, curves turn only for nearly incompressible equations 
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Figure 3. Same as figure 2 but for the results with the N = 1 
equation of state. 



of state, but curves behave monotonically for compressible 
equations of state. Since gravitational radiation carries away 
the angular momentum and the energy, the turning points 
on the curves correspond to critical points where dynamical 
instability sets in. In this sense, our present results differ 
from those of LRS1 and LRS2. 

From the observational point of view, values of the an- 
gular velocity have an important meaning because they are 
related to frequencies of emitted gravitational wave. In par- 
ticular, those at critical points divide the ordered frequencies 
for quasi-periodic stages and the frequencies for dynamically 
inspiraling stages. As seen from figures 1 and 2 the critical 
frequencies obtained from the ellipsoidal approximations are 
several percent larger that exact values of our present irro- 
tational binary models. 
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Figure 4. Physical quantities of irrotational binary sequences 
for several polytropic indices, (a) Total angular momentum as 
a function of a binary separation, (b) Orbital angular velocity 
as a function of the total angular momentum. Different curves 
correspond to different polytropic indices: N = (dash dotted 
line), N = 0.5 (dotted line), N = 0.7 (dashed line) and N = 1 
(solid line). 



These tendencies can be more explicitly shown in figure 
4. In figure 4, stationary sequences along which the mass, 
the entropy and the circulation are the same, i.e. quasi- 
evolutionary sequences, are shown for polytropic indices 
N = 0, 0.5, 0.7 and 1. As seen from figure 4, evolutionary 
sequences with stiffer equations of state have turning points 
but there are no turning points for more compressible equa- 
tions of state. As discussed above, models corresponding to 
turning points are critical ones beyond which dynamical in- 
stability of the binary stars will set in. As figure 4 shows, 
for the binary stars with index N J> 0.7, there are no critical 
points where dynamical instability will set in before contact 
stages. It is important to note that the equation of state 
of neutron star matter can be approximately represented 
by polytropes with index N = 0.5 ~ 1 (see e.g. Shapiro & 
Teukolsky 1983). Thus, as the separation decreases due to 
gravitational wave emission, the inviscid NS-NS binary will 
stably evolve to a contact configuration. 



4 DISCUSSION AND CONCLUSION 

In order to investigate evolutions of binary stars, several crit- 
ical states should be taken into consideration. They are (1) 
a dynamically unstable state, (2) a Darwin-Riemann limit 
where the mass overflow from one component to the other 
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will set in, and (3) a contact phase of binary stars. There can 
be varieties about the final fates of the binary neutron star 
system depending on the order of appearance of these crit- 
ical states as the separation decreases due to gravitational 
wave emission. In particular, for the equal mass binary sys- 
tem, states of the contact phase and of onset of dynamical 
instability are expected to exist. Our results show that for 
the binary system with N <; 0.7 the dynamical instability 
point does never appear on the binary sequence but that 
binary stars do contact before dynamical instability sets in. 
If merging after a contact phase proceeds qua si-statically, 
a dum bbell-li k e configuratio n will be formed ( Eriguchi & 
Hachisx 1985; Hachisu 198C). Since this configuration will 
continue to evolve due to gravitational wave emission, it 
will become dynamically unstable at a certain point on this 
dumbbell sequence and begin to collapse violently because it 
is highly likely that there is a turning point on the dumbbell 
sequence. 

Although our analysis has been made in the framework 
of Newtonian gravity, there is a possibility that the same 
situation will occur even for general relativistic treatments. 
If it will be the case for general relativity, the scenario of 
evolution of irrotational binary systems should be changed. 
In particular, we should consider what will happen after 
contact phases of binary stars. As for non-equal mass binary 
systems, there is a possibility that a similar situation will 
happen, i.e. that the Darwin-Riemann limit or the contact 
phase will occur prior to dynamical instability. 

Consequently, for further quantitative research on real 
binary neutron star systems, it is important to take into ac- 
count the general relativistic effect. General relativity is apt 
to make the radius smaller and the density distribution to 



more contrally condensed. The typical relativistic effect is 



the existence of the innermost stable circular orbit (see e.g. 
Lincoln & Will 1990; Kidder, Will & Wiseman 1992). How- 
ever, as Lai, Rasio & Shapiro (1993b) discussed, the radius 
of this orbit is well inside the radius of the onset of dy- 
namical instability due to the hydrodynamical (tidal) effect 
for most binary neutron star systems as far as the quan- 
tity Rc 2 /GM satisfies the relation Rc 2 /GM £ 5, where R 
and c are the stellar radius and the speed of light, respec- 
tively. This implies that for irrotational binary neutron star 
systems contact phases will appear first as the separation 
decreases. In order to investigate general relativistic effects 
on the fate of binary neutron stars quantitatively, the com- 
petition among the tidal effect to elongate the stellar size 
and the relativistic effect to make it compact will be re- 
quired. Recently the synchronously rotating NS-NS binary 
systems in general relativity have been approximately com- 
puted by several authors (Shibata 1994, 1997; Baumgarte 
et al. 1997). We may be able to extend our method of cal- 
culating irrotational binary configuration in the framework 
of general relativity by using a similar formulation as the 
present Letter. 

The advanced gravitational wave observatory men- 
tioned in Introduction will certainly detect gravitational 
waves from the final phase of NS-NS binary systems in the 
early years of the next century. In order to analyze gravi- 
tational waves and get information from these events, it is 
necessary to compute reliable wave forms at the merging 
stages. Shibata, Nakamura & Oohara (1992) performed dy- 
namical computations of coalescing binary neutron stars and 



reported that the wave form depends on whether the compo- 
nent stars have spins or not (see also Rasio & Shapiro 1996, 
and references therein). Therefore the spin effect plays an 
important role for the investigation of determination of the 
equation of state for neutron star matter by using the wave 
form evaluated from hydrodynamical simulations. However, 
in these dynamical computations they began their computa- 
tions from initial states which were not in exact equilibrium 
of binary stars with spins. In order to have more reliable 
results, it would be desirable to begin dynamical computa- 
tions from exact equilibrium states of spinning binary stars. 
Our new solutions in this Letter can be used as initial states 
for such dynamical computations of coalescing binary stars. 
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